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ABSTRACT 

A  barotropic  and  two-  and  three-layer  baroclinic 
primitive  equation  prediction  models  are  developed  from  a 
ten-layer  model  previously  designed  for  short  range  fore- 
casts in  the  tropics.   Twenty-four  hour  predictions  from 
the  models  are  compared  to  evaluate  any  degradation  in  the 
forecasts  as  the  vertical  resolution  of  the  models  decreases. 
The  barotropic  primitive  equation  model  on  a  global  tropical 
grid  is  used  to  evaluate  the  possibility  of  using  a  numerical 
variational  analysis  as  an  initialization  technique. 
Initialization  of  data  with  the  nonlinear  balance  equation 
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Fig.  1.   Vertical  distribution  of  variables  in  the 
ten-layer  model. 

Fig.  2.   Vertical  distribution  of  variables  in  the 
barotropic  model. 

Fig.  3.   Vertical  distribution  of  variables  in  the 
two -layer  model. 

Fig.  4.   Vertical  distribution  of  variables  in  the 
three-layer  model. 

Fig.  5.   Limited  grid.   There  are  510  computational 
points.   The  grid  spacing  of  1.5  degrees 
longitude  at  22. 5N  is  shown  in  the  upper 
left  corner. 

Fig.  6.    FNWC  Monterey  Tropical  Global  Mercator 

Projection.   There  are  7105  computational 
points  with  grid  spacing  of  2.5  degrees 
longitude  at  the  equator. 

Fig.  7.    Restoration  coefficients  near  the  northern 
boundary  of  the  tropical  global  belt.   A 
persistence  forecast  results  at  the  boundary, 
a  fully  dynamic  forecast  seven  rows  in,  and 
a  blend  in-between. 

Fig.  8.    Twenty-four  hour  250  mb  height  forecast  from 
the  barotropic  model. 

Fig.  9.    Twenty-four  hour  250  mb  height  forecast  from 
the  ten- layer  model. 

Fig.  10.  Total  (T),  zonal  (Z),  and  eddy  (E)  kinetic 
energies  vs.  time  for  the  barotropic  model 
on  the  limited  grid. 

Fig.  11.   Total  (T),  zonal  (Z),  and  eddy  (E)  kinetic 
energies  vs.  time  for  the  two-layer  model 
on  the  limited  grid.   Solid  lines  indicate 
250  mb  data  while  dashed  lines  indicate 
450  mb  data. 

Fig.  12.   Twenty-four  hour  250  mb  height  forecast  from 
the  three -layer  model. 


Fig.  13.  Twenty-four  hour  750  mb  height  forecast  from 
the  three -layer  model. 

Fig.  14.  Twenty-four  hour  750  mb  height  forecast  from 
the  ten-layer  model. 

Fig.  15.   Twenty-four  hour  1000  mb  height  forecast  from 
the  three -layer  model. 

Fig.  16.   Twenty-four  hour  1000  mb  height  forecast  from 
the  ten-layer  model. 

Fig.  17.   Initialization  of  the  500  mb  height  field  by 
the  nonlinear  balance  equation  technique. 
Solid  lines  are  initial  data  while  dashed 
lines  are  the  departure  from  the  5 00  mb 
climatology  used  as  input. 

Fig.  18.   Initialization  of  the  500  mb  height  field 

by  the  numerical  variational  analysis  tech- 
nique.  Solid  lines  are  initial  data  while 
dashed  lines  are  the  departure  from  the  500 
mb  climatology  used  as  input. 

Fig.  19.  Twenty- four  hour  250  mb  height  forecast  from 
the  barotropic  model  on  the  global  grid  with 
restoration  boundary  conditions  and  numerical 

vtiixctoj-uiiaj.    cxii'dj-y  a  j_b     j_iij_  oj-cxi.i.kJ<a,0.LGj.i , 

Fig.  20.   Twenty-four  hour  250  mb  height  forecast  from 
the  barotropic  model  on  the  global  grid  with 
no-flux  boundary  conditions  and  nonlinear 
balance  equation  initialization. 

Fig.  21.   Total  (T),  zonal  (Z),  and  eddy  (E)  kinetic 

energies  vs.  time  for  the  barotropic  model  on 
the  global  grid.   Solid  lines  are  for  the 
basic  model  while  dashed  lines  are  for  the 
model  with  diffusion  terms  included. 

Fig.  22.   Twenty-four  hour  250  mb  height  forecast  from 
the  barotropic  model  with  diffusion  terms 5 
restoration  boundary  conditions  and  numerical 
variational  analysis  initialization. 
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I.   INTRODUCTION 

Weather  prediction  in  the  tropics  has  always  suffered , 
basically  from  a  paucity  of  data.   Hemispheric  numerical 
prediction  models  not  only  suffer  from  this  deficiency,  but 
setting  the  boundary  conditions  in  the  tropics  probably 
also  deteriorates  the  product  in  this  area. 

In  recent  years  there  has  been  some  increase  in 
meteorological  data  over  the  tropical  regions ,  particularly 
in  the  form  of  commercial  aircraft  reports.   This  increase 
allowed  Fleet  Numerical  Weather  Central,  Monterey,  to  commence 
an  objective  analysis  for  a  global  belt  centered  in  the 
tropics   Tv"\  addition  "Fiwn  was  also  investigating  the  use 
of  a  numerical  variational  analysis  scheme  for  this  global 
belt.   Since  this  scheme  simultaneously  solved  for  the 
wind  and  mass  fields  using  the  maximum  amount  of  data,  it 
could  potentially  be  used  as  an  initial  state  for  a  tropical 
numerical  prediction  model,  if  the  fields  were  dynamically 
consistent. 

Multi-layer  baroclinic  numerical  primitive  equation 
prediction  models  for  limited  regions  in  the  tropics  have 
been  developed  by  Krishnamurti  (1969),  Miller  (1969)5 
Harrison  (1969),  and  Elsberry  and  Harrison  (1971).   Vanderman 
and  Collins  (1967)  have  used  a  barotropic  primitive  equation 
model  for  operational  experiments  with  data  from  a  global 
tropical  belt. 
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The  availability  of  tropical  data  for  a  global  belt, 
along  with  the  limited- region  prediction  model  developed  by 
Elsberry  and  Harrison  provided  the  impetus  for  this  work. 
The  object  of  this  thesis  then,  is  to  perform  experiments 
involving  the  initialization  of  data  and  the  development 
of  a  primitive  equation  model  for  short  range  tropical 
forecasts  on  an  operational  basis. 


II.   PRIMITIVE  EQUATION  PREDICTION  MODEL 

A.   BASIC  MODEL 

A  ten- layer  primitive  equation  model,  developed  by 
Harrison  (1969)  and  Elsberry  and  Harrison  (1970,  1971)* 
was  used  to  formulate  less  complex  models  for  possible 
operational  use.   The  original  model  was  developed  for  short 
range  predictions  over  a  3^x15  point  tropical  grid,  with  a 
grid  distance  of  1.5  degrees  longitude  at  22.5  N.   A 
frictionless,  adiabatic  version  of  this  model  required 
approximately  20  minutes  of  computer  time   to  produce  a 
twenty-four  hour  prognosis. 

TV-*,--,     r~,  qi  -,  a  +-  jj  r\-.->  g      ■? 'qv»     -{-Vies     t©Il,l"l£lYSXl    iTlodel     aT*e     ELS     f  <~)"l  1  OWS  ; 

±2  =  -L    (u)   +  fv  -  mlt  +  Fx  (1) 

3  t  d  x 

-It  -  "L  (v)  "  fu  ~  m  ±$-  +  Fy  (2) 

^t  3  y 

-2£   =   -L    (9)    +  HEAT  (3) 

dt 

_La   =   -L    (  q  )    +  MOISTURE  (4.) 

o  ** 

34 1000       T   /x         n  /c\ 

y^ —    =  "L  (*iooo ^  (5) 

|^  =  -m2  [jL  (  a  )  +  a     (  v  j]  (6) 

A<f>   -  @cPA(Io§o-)R^  (7) 

where  2p^        /su\      ^        /svfl         -v 

L(s)  =m  L"iru/+irwJ  +4t(ws)  (8) 
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The  ten-layer  model  averaged  parameters  over  two  points 
only,  in  the  direction  of  the  space  derivative.   Finite 
differencing  on  a  constant  pressure  surface ,    for  the  general 
variable  S  is  as  follows : 

Xj     (SV  )  =  {  [S(I>  J+1;K)1    J   [  V(I,J+1,K)     T(I,J+l,g)1 

*y        m      i    L  J     L    m(JH-l)       m(J)    J 

-fs(I,J-i,K)l   |  [  IILlJ^L]   +  v(l,J-l?K)]|  1 
L  J      L   m(j)         m(J-l)   J  J  Ay 

where 

S(l,J+i,K)  =  S(I,J+1,K)  +  S(I?J?K)  } 

2 

S(I,J-i,K)  =  S(l,J-l,K)  +  S(I,J,K) 

2 

An  important  feature  cf  the  mode]  is  the  vertical  distri- 
bution of  variables ,  as  shown  in  Fig.  1_,  since  vertical 
motion  is  not  averaged.   All  other  parameters  are  averaged 
uniformly  in  the  vertical,  except  at  the  lowest  level. 
Vertical  differencing  is  as  follows : 

2 


3P 


(US)  =  |  [CJ  (I,J,K+1)]  [s(I,J,K+l)] 

-[(J  (I,J,K-1)]  [S(I,J,K-1)]}  1   s 

where 

S(I,J,K+1)  =  S(I,J,K+2)  +  S(I,J,K)  ^ 

2 
S(I,J,K-l)     Sfl,J,K)  +  S(l,J,K-2) 

2 
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Arakawa  (1966)  showed  this  finite  difference  flux  form 
conserves  the  mean  square  of  potential  temperature  and  mean 
square  kinetic  energy  while  preventing  nonlinear  instability 
caused  by  aliasing.   The  model  remained  computationally 
stable  with  minimal  smoothing  or  averaging  utilizing  a  six- 
minute  time  step. 

The  first  time  step  is  forward ,  where  for  the  general 
parameter  S,  time  step,  n,  and  time  difference,  ,At, 

sn=l  =  sn=o  +At   f (sn=o} 

Centered  or  leapfrog  time  steps  are  used  for  subsequent 
time  steps ,  where 

sn+l  =  sn-l  +  2At  f(Sn^ 

Reinitialization,  to  avoid  solution  separation,  during  the 
twehty-±'our  hour  period  was  not  found  to  be  necessary. 

B.   BAROTROPIC  ONE-LAYER  MODEL 

Modification  of  the  ten-layer  model  into  a  one-layer 
barotropic  model  was  the  first  step  in  the  development  of 
an  operational  model.   The  250  mb  surface  was  chosen  as  the 
optimum  forecast  level  due  to  the  relatively  high  density 
of  observations  in  the  tropics  at  this  level,  according  to 
Krishnamurti  (1970). 

Vertical  distribution  of  variables  in  the  one-layer 
barotropic  model  is  shown  in  Fig.  2.   Required  inputs  are 
the  winds  and  geopotential  at  250  mb.   A  ficticious  1000  mb 
geopotential  is  derived,  utilizing  the  barotropic  assumption, 
to  maintain  a  similarity  between  the  barotropic  one-layer 
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and  the  ten- layer  model.   Since  this  is  a  barotropic 

model  Equations  (3)  and  (4)  are  not  required.   In  addition, 

the  frictional  terms  in  Equations  (l)  and  (2)  are  neglected. 

Initial  experiments  were  conducted  using  the  510  point 
grid  that  was  used  in  the  development  of  the  ten- layer 
model.   This  grid  is  illustrated  in  Fig.  5.   Approximately 
two  minutes  of  computer  time  is  required  for  a  twenty-four 
hour  prognosis. 

Later  work  was  performed  using  the  FNWC  (Monterey) 
Tropical  Global  Mercator  Projection,  shown  in  Fig.  6.   This 
is  a  ikkxkg   point  grid  with  a  grid  distance  of  2.5  degrees 
longitude  at  the  equator.   A  twenty-four  hour  barotropic 
prediction  on  the  global  band  requires  approximately  25 

mTrin  +  oo     r»-P    C01TlT)U"he2*    "hnrop 

C.   BAROCLINIC  TWO-LAYER  MODEL 

As  a  second  step,  a  two -layer  model  of  the  atmosphere 
was  developed  to  include  the  effects  of  friction  in  the 
surface  layer.   With  the  exception  of  Equation  (4),  the 
full  set  of  equations  was  used  in  this  model.   However,  in 
addition  to  not  predicting  moisture,  diabatic  heating  was 
neglected  so  that  in  Equation  (3)  only  advection  was  per- 
mitted to  change  the  initial  temperature  field. 

The  vertical  distribution  of  variables  in  the  two-layer 
model  is  shown  in  Fig.  3.   Required  inputs  are  the  winds, 
geopotential,  and  temperatures  at  950  and  250  mb,  and  the 
geopotential  and  temperature  at  1000  mb. 
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The  friction  terms  are  given  by  the  bulk  aerodynamic 
transport  relations 


*x   - 

F      = 

y 

-g 

-g 

)  P 
3p  ' 

where 

3  Cd 

^  Cd 

The  values  of  the  drag  coefficient,  Cd,  are 

Cd  =0.5  |  V|^  x  10  "3, 
for  winds  from  1  to  15  meters  per  second,  and 

Cd  =  206  x  10~3, 
for  wind  speeds  greater  than  15  meters  per  second.  The 
wind  udcJ  in  the  above  calculation  is  a  linear  function 
of  the  950  mb  wind  with  the  direction  held  constant  and 
the  speed  reduced  20  %.  Friction  is  considered  only  in  the 
layer  between  1000  and  900  mb,  where  the  stress  and  the 
vertical  gradient  of  stress  are  both  assumed  to  be  zero 
at  900  mb. 

As  stated  earlier,  the  upper  level  of  most  interest  was 
the  250  mb  level.   However,  this  implies  that  the  wind  field 
in  the  very  deep  upper  layer  is  defined  by  a  level  near  the 
upper  extreme  of  the  layer.   Consequently,  the  two  layer 
baroclinic  model  was  also  tested  in  experiments  with  the 
350,  450,  or  550  mb  level  replacing  250  mb  as  the  upper 
wind  level. 
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All  experiments  with  the  two-layer  model  were  conducted 
using  the  510  point  grid.   A  twenty-four  hour  prognosis 
required  approximately  seven  minutes  of  computer  time. 

D.  BAROCLINIC  THREE-LAYER  MODEL 

Inclusion  of  data  at  the  750  mb  level  produced  a  three- 
layer  version  of  the  baroclinic  model  to  provide  a  more 
realistic  approximation  of  the  atmosphere „   The  800  mb 
above  the  surface  layer  was  divided  into  two  layers  at  the 
500  mb  level.   The  development  of  this  model  also  presented 
the  opportunity  to  evaluate  the  effect  of  increased  vertical 
resolution  versus  the  resulting  increase  in  computer  storage 
and  time  requirements „   Fig.  4  illustrates  the  vertical 
distribution  of  variables  for  the  baroclinic  three -layer 
model. 

All  experiments  with  this  model  were  conducted  with  the 
510  point  grid.   Computer  time  of  approximately  nine  minutes 
is  required  for  a  twenty-four  hour  prediction. 

E.  BOUNDARY  CONDITIONS 

The  computationally  stable  integration  of  the  primitive 
equations  requires  a  set  of  boundary  conditions  that  are 
consistent  with  the  finite  difference  scheme. 

Flux  through  the  top  of  the  domain  is  eliminated  by 
setting   CO  =  0. 
At  the  lowest  levels  1000  mb^  Equation  (5)  is  imposed. 

Experiments  were  conducted  with  three  sets  of  horizontal 
boundary  conditions. 
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1.   No-flux  "boundaries 

Slsberry  and  Harrison  (1971)  discuss  the  boundary 
conditions  which  are  applicable  to  the  basic  model. 
Initialization  was  accomplished  using  the  nonlinear  balance 
equation  after  computing  the  rotational  component  of  the 
wind. 

The  requirement  for  east-west  boundary  conditions 
on  the  510  point  grid  is  eliminated  by  making  the  data 
cyclical  in  the  initialization.   Naturally s    east-west 
boundary  conditions  are  not  necessary  on  the  global  band. 

By  setting 


"jL±i  v,,      u.        =  u,,    and     cj>  =  <f> 


o+i    — —  vj>  -j+i  ~y   ~"~  Tj+i   t  j 


where  .  , -,  is  the  northern -meet  row  ^^   th^  grid,  ^o  flux  of 
mass  or  energy  is  allowed.   The  southern  boundary  is  treated 
in  a  similar  manner. 

20   Restoration  or  constant-flux  boundaries 

Kesel  and  Winninghoff  (1970)  describe  a  constant- 
flux,  restoration  boundary  technique  that  eliminates  both 
the  necessity  for  altering  the  initial  mass  structure  and 
the  false  reflection  of  the  computational  mode  at  the 
boundary.   Initial  height  and  wind  fields  near  the  boundaries 
are  preserved.   At  the  completion  of  each  time  step   the 
newly  predicted  values  of  the  fields  are  restored  toward 
values  at  the  previous  time  step  with  a  specified  field  of 
restoration  coefficients.   Figure  (f)  illustrates  the 
restoration  field  near  the  northern  boundary.   The  procedure 
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results  in  a  fully  dynamic  forecast  seven  rows  in  from  the 
boundary,  a  persistence  forecast  for  the  two  rows  at  the 
boundary s    and  a  blend  in-between „   Restoration  boundary 
conditions  were  restricted  to  experiments  using  the  global 
grid. 

3.   Interpolated  boundary  conditions 

Krishnamurti  (1969)  used  implicit  north  and  south 
poles  to  define  boundary  values.   Geopotential  values  at 
the  north  pole  are   defined  by  the  mean  value  at  row  J3 
where  the  north  boundary  is  row  J+l.   The  three  components 
of  the  wind  are  assumed  to  be  zero  at  the  pole.   A  linear 
interpolation  between  the  pole  and  row  j  is  performed  each 
time  step  to  calculate  the  values  of  the  four  parameters 
at  i.uw  jil.   A  similar  calculation  specifics  the  southern 
boundary.   Experiments  with  interpolated  boundary  conditions 
were  confined  to  the  global  band. 
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III.   DATA  INITIALIZATION 

The  importance  of  initial  conditions  for  a  primitive 
equation  model  has  long  been  recognized.   Raw  meteorological 
data,  used  as  an  input s    quickly  results  in  high  frequency 
inertio-gravity  oscillations.   Shuman  and  Hovermale  (I9^9)j 
Benwell  and  Bretherton  (1968),  Smagorinsky  et  al.  (1965), 
Charney  (1955).?  and  Hinkelmann  (1951)  all  propose  an  initial 
state  with  no  divergence  or  vertical  motion,  to  suppress 
these  oscillations.   Phillips  (i960),  however,  showed  that 
for  a  simple  linear  model  an  initial  divergence  equal  to 
that  suggested  "by  quasi-geostrophic  balance  is  required. 
Hinkelmann  (1959)  demonstrated  that  high  frequency 
oscillations  will  be  filtered  out  if  the  time  derivative 
of  divergence  vanishes  in  the  initial  field,  that  is  if 

1HR  =  0  (9) 

where  D  is  the  divergence  and  n=l,2,35 . . . « . .   Note  that 
using  the  balance  equation  to  specify  the  initial  conditions 
results  in  satisfying  the  above  conditions  for  n=l.   Finally, 
this  condition  has  been  satisfied  for  n=l  and  2  by  Miyakoda 
and  Moyer  (1968)  and  Nitta  and  Hovermale  (1969),  by 
integrating  the  primitive  equations  forward  and  backward 
until  the  initial  field  meets  the  condition.   Two  types  of 
initialization  were  tested  in  this  series  of  experiments. 


18 


Ac   NONLINEAR  BALANCE  EQUATION 

Elsberry  (1969)  utilized  the  nonlinear  balance  equation 
and  special  boundary  conditions  to  formulate  a  diagnostic 
model  that  would  match  the  prognostic  model  (ten-layer  with 
no-flux  boundary  conditions).   The  nondivergent  part  of  the 
wind  was  obtained  from 

where  I  is  the  relative  vorticity  from  the  observed  wind 
components  u  and  v»  Smoothing  the  relative  vorticity  fields 
in  the  vertical  improved  the  vertical  consistency  of  the 
stream- function  fields.  No  normal  component  of  the  wind  was 
allowed  at  the  north  and  south  boundaries  _,  implying  that 
these  boundaries  were  streamlines.  Thus  j n  and  is  were 
set  equal  to  constants  at  each  of  the  boundaries  with 


T„-Ts=*Ay5 


where  u  is  the  mean  zonal  wind  component  averaged  over  the 
grid,  and  A  y  is  the  distance  between  the  north  and  south 
boundaries.  In  the  manner  of  Krishnamurti  (1969)  the  grid 
was  extended  by  five  points  between  the  east  and  west  ends 
of  the  grid  which  are  linearly  matched,  providing  cyclical 
continuity. 

The  Afield  was  used  to  determine  the  geopotential  from 
the  nonlinear  balance  equation } 

T  3  x  ,   9  y 

Boundary  conditions  on  ty  at  the  northern  and  southern  walls 
were  made  geostrophically  consistent  with  the  streamf unction 
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difference  (jn-is)  by  using  a  mean  value  of  f  across  the 
limited  region  grid.   On  the  global  belt  the  averages  of 
the  observed  geopotential  heights  along  the  walls  were  used  for 
the  boundary  conditions. 

B.   NUMERICAL  VARIATIONAL  ANALYSIS 

Sasaki  (1958)  introduced  the  variational  method  to 
obtain  dynamically  consistent  initial  fields,  showing 
theoretically  that  quasi-geostrophically  balanced  or 
divergence-free  fields  could  be  obtained  as  a  solution  to 
the  Euler-Lagrange  equations .   Filtering  high  frequency 
oscillations  using  Equation  (9)  as  a  low  pass  filter,  where 
D  is  not  limited  to  divergence,  but  may  be  any  meteorological 
parameter  to  be  filtered,  was  also  proposed  by  Sasaki  (1970). 
ijOiig  wave  tj^ucL i/.i.i  ma   ci.it  u.£>eu  with  quas  jl-s t;ca-<.y  s  oa"Cc  ccn 
ditions  as  a  dynamical  constraint  in  a  numerical  variational 
analysis.   This  method,  as  adapted  by  Lewis  (1971)  for 
operational  use  by  FNWC,  Monterey,  was  tested  as  an 
initialization  scheme  for  the  barotropic  prediction  experi- 
ments on  the  global  grid. 

The  variational  analysis  was  obtained  by  minimizing  the 
functional 

J  =}f{«[(*-~f  +  (v-rf]  +j  (z-z)2 

where  the  integration  is  over  the  surface  S  being  analyzed. 
Here  u,  v,  and  z  are  the  observations  interpolated  to  grid- 
points  by  an  objective  analysis  scheme,  and  u,  v,  and  z 
are  the  desired  fields  which  result  from  the  minimization 

20 


of  Jo  oC  and  3   are  weighting  functions  dictated  by  the 
density  and  reliability  of  observations.   The  last  two  terms 
in  the  functional  can  be  interpreted  as  a  control  on  the 
amount  of  acceleration  in  the  analysis.,  and  are  related  to 
the  low  pass  filter  characteristic  mentioned  earlier.   The 
larger  the  ratio  of  the  dynamical  weight  °t>   to  the  obser- 
vational  weights.,  «*  and  p  ,  the  more  the  analysis  will  be 
forced  toward  a  steady  state  where 

-lu  =  -  -I  ~  -V-Vu  +  f v  =  0 

at       s  9x 

and 

JLX  =  -  1  5*2  -V.nv  -  fu  =  0. 
at     $    ay 

Steady  state  can  only  be  attained  in  the  limits  if  °/c< 
and  «/^  are  -ero  ove3  the  whole  field.   The  inclusion  of 
some  small  tendencies  in  the  fields  however,  may  be  desirable 
for  the  initialization  of  a  primitive  equation  prediction 
model. 

To  enhance  the  operational  value  of  the  variational 
analysis  scheme,  Lewis  (1971)  developed  a  technique  to 
avoid  solving  the  three  equations  for  the  u,  v,  and  z  fields 
simultaneously.   However,  this  technique  involved  grouping 
the  inertial  terms  with  the  local  change  terms,  and  minimizing 
the  imbalance  between  the  remaining  terms,  the  pressure 
gradient  and  coriolis  terms.   The  usefulness  of  the  variational 
analysis  for  initialization  was  degraded  by  this  technique, 
particularly  in  high  wind  regions .   The  technique  was  then 
modified  to  include  an  approximation  of  the  inertial  terms, 
by  evaluating  them  using  the  input  fields,  u  and  v\ 
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IV.   RESULTS 

A.   PREDICTION  EXPERIMENTS  WITH  THE  LIMITED  GRID 

This  group  of  experiments  is  an  attempt  to  test  how  much 
a  barotropic  model  and  simple  baroclinic  models  degrade  the 
results  of  the  ten-layer  model.   Data  initialization  was 
accomplished  using  the  nonlinear  balance  equation  technique 
and  the  prediction  models  had  no-flux  boundary  conditions. 
The  input  data  contained  both  an  upper  level  cold  low  and 
the  anticyclonic  flow  over  an  incipient  tropical  storm. 
1.   Barotropic  Model 

The  barotropic  model  provided  a  twenty-four  hour 
prognosis  that  compared  very  favorably  with  the  twenty-four 
hour  results  from  the  ten-layer  model.   Figures  8  and  9 
illustrate  the  twenty-four  hour  prognoses  from  the  respec- 
tive models.   The  barotropic  prognosis  shows  essentially 
the  same  height  pattern  as  the  ten-layer  prognosis.   The 
greatest  difference  is  in  the  height  of  the  cold  low, 
where  the  barotropic  prediction  is  approximately  ten  meters 
higher  than  the  ten-layer  prediction. 

Figure  10  illustrates  the  change  of  the  total, 
zonal,  and  eddy  kinetic  energies  with  respect  to  time  for 
the  barotropic  model.   Although  an  apparently  spurious 
increase  in  the  eddy  kinetic  energy  is  noted,  no  degradation 
of  the  twenty-four  hour  prognosis  was  noted. 
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2.   Two-layer  Model 

The  attempt  to  predict  the  250  mb  surface  with  the 
two-layer  model  was  not  successful.   The   exponential 
growth  of  eddy  kinetic  energy  associated  with  this  failure 
is  shown  in  Fig.  11,   where  the  solid  lines  indicate  the 
250  m"b  kinetic  energies.   The  dashed  lines  illustrate  the 
450  mb  kinetic  energies.   Two- layer  models  with  upper  levels 
of  350  and  550  mb  also  did  not  exhibit  unreasonable  growth 
of  the  eddy  kinetic  energy.   This  implies  that  perhaps 
successful  predictions  could  be  carried  out  when  the  broad 
upper  layer  is  defined  with  data  from  levels  lower  than 
250  mb,  or  with  input  data  which  has  been  averaged  over 
several  layers. 

q     nru  — % ~  ^  _  n  oygyi  MO ' 3  6  ] 

Twenty-four  hour  prognoses  from  the  three -layer 
model  were  very  encouraging.   The  twenty-four  hour  250  mb 
prognosis ,  as  shown,  in  Fig.  12,  is  almost  identical  to  the 
ten-layer  prognosis  shown  in  Fig.  9°   There  are  no  appreciable 
differences  in  the  heights „   Of  course,  the  three-layer 
model  also  provides  predictions  at  750,  950,  and  1000  mb. 
The  750  mb  twenty-four  hour  prognoses  as  illustrated  in 
Fig.  13  and  14-  also  have  no  significant  differences.   Fig.  15* 
the  three-layer  1000  mb  twenty-four  hour  prognosis,  while 
having  the  same  pattern  as  Fig.  16,  the  ten-layer  version, 
has  filled  the  tropical  cyclone  approximately  ten  meters 
more,  thus  reducing  the  gradient  surrounding  the  storm. 
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Thus  a  three-layer  model  which  would,  incorporate 
both  the  boundary  layer  and  the  upper  tropospheric  flow,, 
which  are  known  to  be  important  in  the  tropics ^  may  be 
feasible.   The  middle  layer  would  be  a  useful  buffer  and 
add  an  additional  degree  of  freedom  for  baroclinic  develop- 
ment.  The  scarcity  of  data  in  the  tropics  for  this  mid- 
layer  could  lead  to  the  use  of  low  cloud  motions ^  as 
measured  from  the  Applied  Technology  Satellite ,    for 
initialization  of  this  level. 

Unfortunate ly_,  only  one  data  set  was  available  for 
the  experiments  with  the  limited  grid  and  further  experi- 
ments with  more  data  sets  are  required  before  these  results 
may  be  considered  conclusive. 

B.   EXPERIMENTS  ON  THE  GLOBAL  BAND 

Use  of  the  nonlinear  balance  equation  for  the 
initialization  of  the  primitive  equation  model  has  been 
demonstrated  previously  by  Harrison  (1969)  and  Elsberry  and 
Harrison  (1971).   Krishnamurti  (1969)  also  used  this  tech- 
nique for  the  initialization  of  a  primitive  equation  model. 
The  predictions  obtained  from  the  barotropic  and  two-  and 
three-layer  models  tend  to  confirm  these  results.   Although 
only  the  rotational  component  of  the  wind  and  a  smooth  height 
field  from  the  balance  equation  are  used  as  initial  data5 
a  realistic  prediction  of  systems  is  the  resultant . 

It  was  considered  sufficient  to  test  the  initialization 
properties  of  the  numerical  variational  analysis  with  the 
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barotropic  model  in  order  to  save  computer  time.   To  assure 
a  good  data  set,  climatological  fields  were  used  rather 
than  actual  data  in  the  global  band  experiments.   Both  500  mb 
and  300  mb  climatological  fields  were  used.   Although  only 
the  500  mb  results  are  discussed  and  illustrated,  the  same 
type  problems  were  encountered  with  both  fields  and  the 
solutions  were  similar.   Since  the  size  of  the  FNWC  tropical 
Mercator  projection  would  not  permit  a  detailed  view  of  the 
entire  grid  in  this  publication,  only  a  small  portion  of 
the  grid  is  presented.   This  portion  was  chosen  as  being 
representative  of  the  problems  and  solutions  encountered 
with  the  whole  grid  and  includes  the  entire  north-south 
extent  of  the  grid,  the  western  Pacific  Ocean,  the  Asian 
coastline,  Japan,  Indonesia  and  Australia.   The  500  mb 
climatological  field  for  January  in  this  area  features  a 
broad  trough  in  the  north,  a  band  of  maximum  winds  between 
30  and  40N,  flat  high  pressure  from  2 ON  to  20S  and  basically 
zonal  flow  in  the  southern  hemisphere . 

Figure  17  represents  the  nonlinear  balance  equation 
initialization.   Solid  lines  are  the  initial  data  while  the 
dashed  lines  are  the  departure  of  the  initial  data  from 
the  500  mb  climatology.   Smoothing  near  the  boundaries  is 
apparent  as  is  a  large  overall  height  rise  caused  by  the 
initialization  boundary  conditions .   On  a  limited  tropical 
grid,  setting  the  heights  on  the  north  and  south  boundaries 
equal  to  the  mean  height  is  reasonable,  however,  on  the 
global  band,  with  boundaries  well  into  mid-latitudes,  this 
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is  no  longer  true.   In  the  region  of  a  trough  along  the 
boundary,  the  field  will  be  relaxed  toward  a  boundary  that 
is  too  high 5  while  the  opposite  effect  will  occur  in  the 
area  where  a  ridge  is  along  the  boundary.   New  boundary 
conditions  will  be  required  if  the  nonlinear  balance 
equation  initialization  technique  is  to  be  used  on  a  grid 
with  boundaries  in  mid -latitudes . 

Initial  heights  as  determined  by  the  numerical  variational 
analysis  technique  are  shown  in  Fig.  l8  as  solid  lines.   The 
dashed  lines  depict  the  deviation  of  the  initial  heights 
from  the  climatological  500  mb  data.   It  is  readily  apparent 
that  this  technique  does  not  significantly  alter  the  input 
field. 

Perhaps  the  greatest  advantage  to  be  realized  by  using 
the  numerical  variational  analysis  technique  for  initialization 
is  provided  by  the  opportunity  to  vary  the  weighting  functions 
over  the  analysis 0   In  sparse  data  regions  the  dynamical 
weight,  o<  ,  can  be  increased.   In  mid-latitudes  the  height 
field  weight,  f%   ,  can  be  increased,  while  in  the  tropics 
the  wind  field  weight,  ex  }    can  be  increased.   The  proper 
utilization  of  the  weighting  functions  should  guarantee  a 
dynamically  consistent  input  for  primitive  equation  models. 

To  determine  if  the  variational  analysis  technique 
would  initialize  the  barotropic  model,  a  twenty-four  hour 
prediction  with  initialization  by  this  scheme  was  compared 
to  a  similar  prediction  with  initialization  by  the  nonlinear 
balance  equation  technique.   Restoration  boundary  conditions 


26 


were  used  in  the  model  with  the  variational  initialization. 
Figure  19  illustrates  the  twenty-four  hour  prediction  for 
this  technique.   The  twenty-four  hour  prediction  for  the 
model  initialized  with  the  "balance  equation  is  shown  in 
Fig.  20.   This  model  used  no-flux  boundary  conditions. 
Ignoring  the  height  differences  caused  by  the  initialization 
techniques,  and  particularly  in  regions  away  from  the 
boundaries ,  the  predictions  are  very  similar.   Only  the 
developing  short  wave  features  in  the  wind  maximum  between 
30N  and  40N  are  disturbing.   In  all,  three  sets  of  boundary 
conditions ,  no-flux,  restoration,  and  interpolated,  were 
tested  with  the  two  initialization  techniques .   All  combi- 
nations developed  the  short  wave  features  mentioned  above. 
There  ~^ztz   no  problems  in  the  Dredoiiiinaiitly  zonaJ  southern 
hemisphere,  nor  were  there  any  short  period  height  oscillations 
at  the  equator  as  noted  by  Winninghoff  (1971)  and  Krish- 
namurti  (personal  communication).   The  conclusion  to  be 
drawn  from  these  limited  experiments  then  is  that  the  short 
wave  features  are  independent  of  either  boundary  conditions 
or  initialization  and  that  the  variational  analysis  technique 
will  successfully  initialize  the  barotropic  primitive 
equation  prediction  model 0 

The  solid  curves  in  Fig0  21  illustrate  kinetic  energy 
changes  for  the  restoration-variational  analysis  prediction 
model  for  the  forecast  period.   The  growth  of  kinetic  energy 
is  attributed  to  an  accumulation  of  energy  in  the  smallest 
permitted  scales  as  the  eddy  kinetic  energy  appears  to  be 
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the  forcing  function.   This  phenomena  was  also  noted  by 
Kesel  (1968),    in  a  barotropic  primitive  equation  experiment 
using  a  hemispheric  grid. 

The  accumulation  of  the  eddy  kinetic  energy  was  controlled 
by  adding  a  term  to  the  momentum  equations  to  simulate 
lateral  diffusion  of  momentum,  in  the  form 

yC   m  \7  u  and  Km  \7   v* 

where  )X.  ,    the  lateral  eddy  diffusion  coefficient  is  a 
function  of  the  grid  distance. 

A  twenty- four  hour  prognosis  as  obtained  from  the 
barotropic  model,  with  lateral  diffusion  terms  included 
(K=5  x  10  )>    restoration  boundary  conditions  and 
variational  analysis  initialization  is  shown  in  Fig,  22. 
Note  that  the  aforementioned  short-wave  features  are  absent. 
The  dashed  lines  in  Fig.  21  illustrate  the  kinetic  energy 
changes  for  the  model  with  the  diffusion  term  included. 
That  the  lateral  diffusion  term  prevents  accumulation  of 
energy  in  the  small  grid  scales  is  evident. 
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V.   CONCLUSIONS  AND  RECOMMENDATIONS 

A  barotropic  or  three-layer  primitive  equation  prediction 
model  for  use  in  the  tropics  on  an  operational  basis  seems 
well  within  reach,  whether  on  a  limited  region  or  a  global 
basis.   The  basic  model  for  use  on  the  FNWC  tropical 
Mercator  grid  has  been  tested.   Further  experiments  are 
required  to  isolate  the  cause  of  and  either  eliminate  or 
control  the  accumulation  of  eddy  kinetic  energy.   As  the 
prediction  models  become  more  refined,  a  staggered  grid 
system  should  be  investigated  to  further  reduce  computer 
time  and  storage  requirements. 

The  successful  use  of  the  nonlinear  balance  equation  for 
initialization  of  a  limited-region  primitive  equation 
prediction  model  has  been  demonstrated  by  Krishnamurti  (19^9) 
and  Elsberry  (1969).   However,  if  this  scheme  is  to  be 
utilized  when  the  boundaries  are  in  areas  with  height  gra- 
dients larger  than  those  encountered  in  the  tropics,  more 
suitable  boundary  conditions  must  be  developed.   Over  a 
global  band,  no-flux  boundary  conditions  may  not  be  crucial 
and  other  boundary  conditions  may  be  applicable. 

The  numerical  variational  analysis  technique  has  been 
adapted  for  use  in  initializing  the  barotropic  primitive 
equation  model.   Using  the  nonlinear  balance  equation 
initialization  scheme  as  a  comparison,  it  would  appear  as 
if  the  numerical  variational  analysis  will  successfully 
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initialize  a  primitive  equation  prediction  model „   Further 
experiments  with  operational  data  will  be  required  to 
determine  the  optimum  values  of  the  weighting  functions . 
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Fig.  1.   Vertical  distribution  of  variables 
in  the  ten- layer  model. 
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Fig.    2„      Vertical  distribution   of  variables 
in  the   barotropic  model. 
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Fig,  3.   Vertical  distribution  of  variables 
in  the  two- layer  model. 
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Fig.  4.   Vertical  distribution  of  variables 
in  the  three -layer  model. 
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Fig.  5.   Limited  grid.   There  are  510  computational 
points..   The  &y\  c\   spacing  of  i  S  degrees 

IGxj^jl  ouut     ao    22,5™     J-t5     hiJU'.*xi     j.j  t     uie     upper 

left   corner. 
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Fig.  6.   FNWC,  Monterey,  Tropical  Global  Mercator 
Projection.   There  are  7105  computational 
points  with  grid  spacing  of  2.5  degrees 
longitude  at  the  equator. 
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Fig„  7.   Restoration  coefficients  near  the 

northern  boundary  of  the  tropical  global 
belt.   A  persistence  forecast  results  at 
the  boundary,  a  fully  dynamic  forecast 
seven  rows  in,  and  a  blend  in-between. 
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Fig.  10.  Total  (T),  zonal  (z),  and  eddy  (E)  kinetic 
energies  vs.  time  for  the  barotropic  model 
on  the  limited  grid0 
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Fig.  11.   Total  (T)3  zonal  (z),  and  eddy  (E) 

kinetic  energies  vs .  time  for  the  two- 
layer  model  on  the  limited  grid.   Solid 
lines  indicate  250  mb  data  while  dashed 
lines  indicate  450  mb  data. 
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Fig.  17.   Initialization  of  the  500  mb  height  field  by  th« 
nonlinear  balance  equation  technique.   Solid 
lines  are  initial  data  while  dashed  lines  are 
the  departure  from  the  500  mb  climatology  used 
as  input . 
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Pig.  18 „   Initialization  of  the  500  mt>  height  field  "by  the 
numerical  variational  analysis  technique.   Solid 
lines  are  initial  data  while  dashed  lines  are  the 
departure  from  the  5 00  mb  climatology  used  as  input 
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Fig.  19.   Twenty-four  hour  250  mb  height  forecast  from  the 

barotropic  model  on  the  global  grid  with  restoration 
boundary  conditions  and  numerical  variational 
analysis  initialization. 
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Fig.  20.   Twenty- four  hour  250  mb  height  forecast  from  the 
barotropic  model  on  the  global  grid  with  no-flux 
boundary  conditions  and  nonlinear  balance 
equation  initialization. 


47 


150- 


T 


125- 


100- 


75- 


50- 


0 


0 


6 


12 
HOURS 


^8 


24 


Fig.  21.   Total  (T),  zonal  (z),  and  eddy  (E)  kinetic 

energies  vs.  time  for  the  barotropic  model  on  the 
global  grid.   Solid  lines  are  for  the  basic  model 
while  dashed  lines  are  for  the  model  with  diffusion 
terms  included. 
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Fig.  22.   Tv/enty-four  hour  250  rab  height  forecast  from  the 

barotropic  model  with  diffusion  terms 3    restoration 
boundary  conditions  and  numerical  variational 
analysis  initialization. 
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